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ABSTRACT 

By performing two-dimensional special relativistic (SR) magnetohydrody- 
namic simulations, we study possible signatures of gravitational waves (GWs) 
in the context of the collapsar model for long-duration gamma-ray bursts. In our 
SR simulations, the central black hole is treated as an absorbing boundary. By 
doing so, we focus on the GWs generated by asphericities in neutrino emission 
and matter motions in the vicinity of the hyperaccreting disks. We compute nine 
models by adding initial angular momenta and magnetic fields parametrically to 
a precollapse core of a 35M Q progenitor star. As for the microphysics, a realistic 
equation of state is employed and the neutrino cooling is taken into account via 
a multiflavor neutrino leakage scheme. To accurately estimate GWs produced 
by anisotropic neutrino emission, we perform a ray-tracing analysis in general 
relativity by a post-processing procedure. By employing a stress formula that in- 
cludes contributions both from magnetic fields and special relativistic corrections, 
we study also the effects of magnetic fields on the gravitational waveforms. We 
find that the GW amplitudes from anisotropic neutrino emission show a mono- 
tonic increase with time, whose amplitudes are much larger than those from 
matter motions of the accreting material. We show that the increasing trend of 
the neutrino GWs stems from the excess of neutrino emission in the direction 
near parallel to the spin axis illuminated from the hyperaccreting disks. We point 
out that a recently proposed future space-based interferometer like Fabry-Perot 
type DECIGO would permit the detection of these GW signals within m 100 
Mpc. 
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Subject headings: massive stars: collapse — gamma-ray bursts — gravitational 
waves — neutrinos — hydrodynamics 



Introduction 



Gamma-ray bursts (GRBs) are one of the most energetic phenomena in the universe. 
Thanks to Swift observation^, it has now become evident that GRBs are basica l ly cat - 



egorized into two, namely short-hard and long-soft bursts (e.g., iHjorth Sz Blooml ( 1201 ll ); 
Nakarl (120071 ) for recent reviews) . More surprisingly, GRBs with some mixe d featu res of 
the two types have been reported (e.g., iGehrels et al.l ( 120061 ). iGal-Yam et al.l (120061 )) pos- 
sibly necessitating a new classification ( iLii et al.l 120101 ) . The mystery of their central en- 
gines seems to be thickening, which has l ong p uzzled astrophysicists since the accidental 
discovery in the late 1960s (see iMeszarosi (120061 ) for review). Regarding the long-duration 



GRB s (LGRBs), ro b ust a ss ociations of the und e rlying supernovae w i th a handful of L GRB s 
e.g.. iGalama et al.l ffl998h : iHiorthet all J2003h : Istanek et al.l J2OO3I): iMalesani et al.l Jiool ; 



Modiaz et a. 



d2006h :l Pian et al. and collective references in lWoosley fc Blooml (120061 ) 



Zhangl (1201 lh 1 and the fact tha t their host galaxies are typically irregular with intense star 



formation ( iFruchter et al.ll2006l ) suggest that they are likely related to the deaths of massive 



(Wooslev 


1993; 


Paczvnski 


1998: 



In the scenario, the collapsed iron core of a massive star forms a tempora ry disk around 



a few Mp) black hole (BH) and accretes at a high r ate (~ 0.1 — 10M^/s, e.g.. IPopham et al. 
dl999[ ): lDi Matteo et alJ(bo02h : lKohri et al.l d2005h : lchen & Beloborodovl (ho07f ): IZalamea & Beloborodov 



(120 111 ) and refer ences there i n) , wh ose g ravitational binding ener gy is the driving source of the 



central engine. iPaczynskil (Il990l ) and iMeszaros fc Reed (119921 ) pioneeringly proposed that 



pairs of neutrino and anti-neutrino illuminated from the hyperaccreting disks that annihilate 
into electron and positron (e.g., v + v — > e~ + e + , hereafter "neutrino pair annihilation") 
can supply sufficient energy to launch GRB outflows by heating material in the polar funnel 
regions. In addition, it is suggested that the strong magnetic fields in the cores of order 
of 10 15 G play also an active role both for driving the magneto-driven iets and for extract- 



ing a significant amount of energy 


Torn the central engine (e.g., ] 


31andford & Znaiek (1977h: 


Mizuno et al. (2004a b): McKinnev 


(2006 


): McKinnev & Naravar 


(2007a b): Komissarov & Barkov 


(2007 


); 


Komissarov & McKinnev (z 


!007); 


Komissarov & Barkov ( 


2009): Barkov & Komissarov 


(2008 


); 


Nagataki (2009h. and references therein). 



http://www.swift.psu.edu/ 
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Al though various possibilities including magnetar models (e. g. , iDai fc Lul (1l998h ; iThompson et al 



( 120041 ) i lUzdensky fc MacFadyenl (120071 ) ; iBucciantini et al.l (120071 )) have been proposed so far, 
there has been no direct evidence to pin down the mechanism of the central engine. This 
is mainly because it is difficult to extract the information from conventional astronomy by 
electromagnetic waves, since high-energy photons are absorbed through interactions in the 
source and by the photon backgrounds. Alternatively, gravitational waves (GWs) are ex- 
pected to be a primary observable to decipher the mechanism of the engine, because they 
imprint a live information hidden deep inside the stellar core and they carry the information 
directly to us without being affected in propagating to the earth. 

Currently long-baseline laser interferometers such as LIGO (lAbbott et al. 1 12005^ VIRGcR 
GE O60C§, and TAMA3 00 (Undo fc the TAMA Collaboration! hoosh are operational (see, 



e.g., iHough et al.l ( 120051 ) for a recent review). For these detectors, core-collapse supernovae 
(CCSNe) have been proposed as one of the most plausible GW sources, therefore an ex- 
tensive study of the GW prediction s ba s ed on sophisticated numer ic al mode l ing h as been 
carried out so far (see, for example, lOttl (120091 ); iFryer &: New! (1201 ll ); iKotakd (120111 ) for re- 
cent reviews). It is noted however that most of them have paid attention to CCSNe that 
leave behind neutron stars (NSs) after explosions. For a reliable prediction of GWs from 
CCSNe, one needs to perform multi-D hydrodynamic simulations equipped with a precise 
neutrino transport scheme which follows the dynamics starting from stellar core-collapse, 
core-bounce, through shock-stall and subsequent growth of hydrodynamic instabilities, the 
neutrino-driven shock revival, to stellar explosion in a consiste nt manner. This is one of the 
most challenging subjects in computational astrophysics (e.g. 



Janka et al.l (120071 )). 



But the numerical modeling to test the collapsar scenario could be much more demand- 
ing. One needs to trace a new path that bifurcates from the above story after bounce, namely 
to the BH formation (phase 1), evolution of the surrounding accretion disk including energy 
deposition to the polar funnel region by neutrinos and/or magnetic fields (phase 2), to the 
launching of the fireballs (phase 3jj. This apparently necessitates the multidimensional mag- 
netohydrodynamic (MHD) simulations not only with general relativity (GR) for handling the 
BH formation, but also with the multi-angle neutrino transfer for treating highly anisotropic 
neutrino radiation from the disks. In the business of CCSN simulations, the most up-to-date 
simulation^ can now follow the multi-angle neutrino transport but limited to a Newtonian 



2 http://www. ego- gw.it/ 
3 http://gco600. aei.mpg.de/ 

4 Here for convenience we call each stage as phase 1, 2, etc. 
5 assisted by accelerating computer powers 
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case JHubenv fe Burrowsl I2OO7L btt et alJl20o3 ISumiyoshi fe Yamadall2012h . or handle GR 



with a sophisticated neutrino transport f lMiiller et al.ll2010l ) but not applicable to a rapidly 
rotating case (due to the assumption of conformal flatness). 

Various approximate approaches have been therefore u ndertaken in the b usiness of 
the collapsar simulations. In the phase 1, GR simulations (IShibata et al.l 120061) upd ated 
to implement a neutrino cooling have reported recently (ISekiguchi &: Shibatal 1201 ll ). in 
which the dynamics after the BH formation to the formation of accretion disk was first 
consistently followed. The numerical studies of the phase 2 are concerned with the sub- 
sequent evolution of accretion disk and the outflow formation in the polar funnel region 
till the jets become mildly relativistic. The centr al BH has been tradit i onally treated by 



a fixe d metric technique in GR simulations (e.g.jMizuno et al.l (12004a); iDe Villiers et al. 



( 20051) : Hawlev fc Kroliki (120061): McKinnev fc Naravan 



ulations (e.g. 



Nagataki et al. 



j2007bh 



Komissarov fc McKinnev 



J2007h : iBarkov & Komissarov! d2008h\ o r by a n absorbing b o unda ry in the Newtoni a n sim 



Harikae et al. 



MacFadyen fc Wooslevl (Il999h: iProga et all ((20031); iFuiimoto et al.l ((20061); 
(Lopez-Camara et al.l ((20091 12010| ) or special relativistic si mulations 



( 20091) ) . As for the microphysics i n these simula t ions, except for lFujimoto et al 



(120061 ) ; iNagataki et al.l (120071 ) ; lHarikae et al.l fl2009|) ; IShibata et al.l (I2007D , a realistic nuclear 
equation of state (EOS) has been replaced by a very phenomenological one (like a gamma-law 
or polytrope) and the neutrino cooling (and heating) has been often neglected for simplicity. 
Numerical studies of the phase 3 are mainly concerned with the dynamics later on, namely, 



the jet propagation to the brea 
polar f unnel region (see, e.g 



tout from the star , by assuming a man u al energy input to the 



Alov et all ( 120001 ): IZhang et al.l ( 120031 ): lLazzati fc Begelman 



j2009h : lNagakura et al.l fcoilal lbh and references therein) All of the studies mentioned ab 



ovc 



may be regarded as complimentary in the sense that the different epochs are focused on, with 
the different initial conditions for the numerical modeling being undertaken. 



To the best of our knowledge, lOtt et al.l ((20111 ) is the only work which extracted the 
GW signals based on their collapsar simulations (in the phase 1). Based on their three- 
dimensional (3D) GR simulations of a 75 M Q star with the use of a polytropic EOS, they 
pointed out that the significant GW emission is associated at the moment of the black hole 
formation, which can be a promising target of the advanced LIGO for a Galactic source. 
In contrast to such a paucity of the GW predictions based on numerical simulations of col- 
lapsars, a number of semi-analytical estimates have been reported so far, which predict a 



signifi cantly strong GW emission due to possible density- inhomogeneities dMineshige et al. 



2002h 

(3); 

d2007h : 



jar or fragmentation instabilitie s in th e collapsar's accretion torii (e. g. 



Davies et al. 



j2002l ) : |Frver et al.l J2OO2I ): iKobavashi fc Meszarosl J2OO3I) 



Corsi k Meszarosl (l2009[ ). and c ollective references inlFrver fc Newl (j201ll )). and the 



precession of the disks due to GR effects ((Romero et al.ll2010l : ISun et al 



van Putten 



Piro fc Pfahl 



2012 1 The predicted 
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GW amplitudes are typically high enough be visible to advanced-LIGO class detectors for a 
100 Mpc distance scale, which is about four order s-of-magnitudes larger than the numerical 
estimate at the black formati on flOtt et al.l 12011). In addition to these GWs produced by 
non-spherical matter motions, Hiramatsu et al.l ( 2005 ); Suwa fc Murase ( 2009 ) pointed out 
that anisotropic neutrino emission from accretion disk could be the source of GWs from 
collapsars, which was originall y proposed as an equally important GW source to the matter 
GW in the context of CCSNe (jEpsteinlll978l ). Since these GWs from collapsars would be a 
smoking-gun signature of the central engine in coincident with the conventional electromag- 
netic messengers as well as neutrino^], it will be very important to put forward theoretical 
predictions of GW signals based on the collapsar simulations, as has been done in the business 
of CCSNe. 

In this work, we study possible GW signatures in the hyperaccreting collapsar disks^ 
by performing two-dimensional special relativistic (SR) MHD simulations of accretion torii 
around a black hole. We compute nine models by adding angular momenta and magneti c 
fields parametrically to a precollapse core of a 35M progenitor star (jWoosley fc Hegerll2006l ). 
As for the microphysics, a realistic equation of state is employed and the neutrino cooling 
is taken into account via a multiflavor neutrino leakage scheme. In our SR simulations, the 
central black hole is treated as an absorbing boundary. By doing so, we focus on the GWs 
generated by asphericities in neutrino emission and matter motions in the vicinity of the 
accretion disks. To accurately estimate GWs produced by anisotropic neutrino emission, 
we perform a ray-tracing analysis in GR by a post-processing procedure. By employing a 
stress formula that includes contributions both from magnetic fields and special relativistic 
corrections, we study also the effects of magnetic fields on the gravitational waveforms. Then 
we discuss their detectability by performing a spectrum analysis. 

The paper opens up with descriptions of the initial models and numerical methods 
(section The main results are given in Section [31 We summarize our results and discuss 
their implications in Section HJ 



6 Especially for a nearby GRB source, e.g.. lAndo et al.l (|2005l ). 
7 which corresponds to evolution in the phase 2 
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2. Numerical Methods and Models 



2.1. Initial Models 




Fig. 1. — Profiles of the specific angular momentum for model 350C (jWoosley &: Heger 
(120061 ). labeled by J350C and the one amplified by a factor of 2 (labeled by 2j 35 oc) and for 
models JO. 6, JO. 8, J1.0, and, J1.2 (from bottom to top), respectively. 



Regarding our precollapse model, we take a 35 star (model 350C) in lWoosley fc Heger 



(120061 ) that is supposed to be one of the most promising GRB progenitor models. To study 
the effects of rotation systematically, we take the following rotational profiles as 

n x* + a n lso (M(x))x 4 



1) 



X 4 + X 4 

where a, f2 , and X are model parameters. M(X) is the mass coordinate at the cylindrical 
radius (X = rsin^), and f2i so is given by Q\ so = j\ so /X 2 , where 7i so is the specific a ngula r 



mome n tum in the last stable orbi t of th e Schwarzshild BH (e.g., iBardeen et al.l ( 119721 ); 



Progal (120051 ); iLopez-Camara et al.l (120091 )). By changing a in the range of 0.6 < a < 1.2, 
we compute 7 non-magnetized models of a = 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, and 1.2, by which 
they are labeled as models JO. 6, JO. 7, and so on (see Figure 1). To see the effects on magnetic 
fields, we compute two more models, in which the initial magnetic field (Bq = 10 10 G or 10 11 
G) are added to model JO. 8 (model J0.8B10 or model J0.8B11), in which the initial field is 
assumed to be purely poloidal and also assumed to be uniform and parallel to the spin axis. 
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Figure 1 shows profiles of the initial angular momentum for some models. As shown, 
our models are taken to possess much more rapid rotation especially in the range of 1 to 4 
M Q in the mass coordinate compared to the original profile of model 350C (black dashed 
line). We set such a rotational profile otherwise the accretion disk cannot survive later than 
~2 s after the onset of gravitational collapse in our simulation. Later on, the accretion disk 
is swallowed to the central object due to the neutrino cooling which deprives the pressure 
support in the disk. In this case, one cannot account for the duration of long bursts, the 
activity of which is typically longer than ~ 2 s, and can last up to tens of minutes. So we 
experimentally choose to adjust the initial angular momentum as in Equation (1). 

Although our models do have higher initial angular momentum, the deviation may not be 
so serious (compare the light blue line that is for the angular momentum amplified by a factor 
of 2 in the original progenitor) considering uncertainties in stellar evolution calculations (such 
as in the treatment of mass- loss, weak interactions, and fluid instabilities (e.g., convection, 
semi convection, rotation, and magnetic fields)). To mimic the original progenitor structure, 
we take Xq to be the size of the Fe core (~ 3000 km). By setting a = 0.8, Qq ~ 1 rad s" 1 , 
the above profile becomes most close to the original profile for the mass coordinates larger 
than 5M Q . As a side- remark, a most prevailing way is to tune the initial angular momentum 
to satisfy its local c entrifugal force to be several percent-levels of the local gravitational 



binding energy (e.g., iMacFadyen fe Woosleyl (Il999h ) . Our modeling may be more realistic 



in the sense that the assumed initial angular momentum profiles capture the basic trend 
obtained in the current progenitor models. 



2.1.1. Hydrodynamics 

As already mentioned in the introduction, a number of Newtonian and SR collapsar 
simulations in the phase 2 have conventionally focused on the evolution of the collapsar disks 
by treating the central BH as an absorbing boundary. Following the tradition of setting the 
boundary from the beginning of the simulation, we also initially impose an absorbing inner- 
boundary condition at the radius of max (10 km, 2r g ) with r g is the Schwarzschild radius 
that is estimated by the accumulating mass inside the inner-boundary. Then we solve the 
dynamics outside the inner boundary up to the outer boundary of the computational domain 
(300 00 km in radiu s ) by o ur SRMHD code assuming axisymmetry and equatorial symmetry 



sec 



Harikae et al.l ( 120091 ) for more details). In our 2D simulations, spherical coordinates 
are employed with logarithmic zoning in the radial direction (r) and regular zoning in the 
polar direction (0 < 9 < 7r/2). The computational domain is covered by 300( r) x 40(8) 



mesh points. Regarding the microphysics, a realistic nuclear equation of state by IShen et al. 
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( 11998 ) is included and t he ne u trino cooling; is treated by a multi-flavor lea k age sch eme 



Epstein fc Pethickl (11981( 1: iRosswog fc Liebendorferl ((20031); iKotake et all (I2003af ). see 
Takiwaki et al.l ( 120091 ) for more details). The gravitational potential is estimated by the 



sum of the Paczynski & Witta-type potential which mimics the gravitational pull from the 
central BH and the self-gravity of material outside t he excised reg i on th at is determined 
by the Poisson equation (see equations (5) and (6) in lHarikae et al.l (120091 ) ). In this paper, 
we examine numerical model s without viscos i ty lik e the alpha prescription (see , however, 



MacFadven fc Wooslev! (119991 ): iLindner et al.l (I2010h : lLopez-Camara et all (120101 )1 



2.2. Extraction of Gravitational Waveforms 



2.2.1. GWs from matter motions 



To extract the gravitational wav eforms from matter motio n s and magnet i c field s, we 
employ the stress formulae derived in iTakiwaki &: Kotakd ( 120111 ) ; IKotake et al.l (120041 ) . For 
convenience, we shortly summarize them in the following. 

In our axisymmetric case, the n on- vanishing quadrupole term is only the plus mode 
(h + ) in the metric perturbation (e.g., iMonchmeyer et al.l ( 119911 )). which can be written as, 



(2) 



Thornd ( 119801 )). The quadrupole amplitude of matter GWs; ^4fo matter that consists of 



(e.g. 

the following three parts, 



A E2 

^20 (matter) 



A E2 , aE2 , aE2 

^20 (hyd) + ^20 (grav) + ^20 (mag) ' 



(3) 



The first term in Equation ([3]) represents the contribution from non-spherical hydrodynamic 
motions, which is expressed by 

G32vr 3 / 2 el 



A E2 
A 20 (hyd) 



da I r dr 
V15 Jo Jo 

xp,W 2 (v r 2 (3fi 2 -l)+v e 2 (2- 3/i 2 ) - V - 6v r v e ^l-u 2 ), 



in which is the effective density defined as, 

e + p+ \b\ 2 



P + 



(4) 



(5) 



Here p, e, p, c, and G, denotes the baryon density, internal energy, pressure, the speed of 
light, and gravitational constant, respectively, and \b\ 2 = b^b^ is related to the energy density 
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of the magnetic fields with representing the magnetic field in the laboratory frame (e.g. 



Takiwaki et al.l (j2009h ). W = — v k v k is the Lorentz boost factor with v k denoting the 

spatial velocity in the spherical coordinates (i — r, 9,(p). \i — cos 9 is a directional cosine. 
The second term in Equation Q represents the contribution from the gravity as, 

G32tt 3 / 2 fX 



A E2 

^20 (grav) 



X 



dp I r dr 
Jo 



ph(W 2 + (v k /c) 2 ) + - 
c z 



X 




((by + ib k f) 



(6) 



-rd r $(3p 2 - 1) + 3<%$ fiy/1 - p? 

where $ denotes the gravitational potential of the self-gravity. The last term in Equation 
([3]) is the contribution from the magnetic fields as, 

G32tt 3 / 2 rl 



A E2 
^20 (mag) 



c 4 V15 

2/o .,2 



^0 



r dr 



x [6 r 2 (3/i 2 - 1) + b e 2 {2 - 3fi 2 ) - b/ - 6b r b ef i^l - /j?]. 



(7) 



Here, we write the total gravitational amplitude as follows for later convenience, 



"-(matter) 



n (hyd) 



(mag) 



""(grav) ) 



where the quantities of the right hand of the equation are defined by combining Equations 
02]) and (|3]) with Equations @, (|6]), and ([7]). Note that by dropping 0{v/c) terms, the 
above formulae re duce to the conven ti onal quadrupole formula employed in the Newtonian 



simulations (e.g., iMonchmeyer et al.l ( 119911 )). In the collapsar disk that we pay attention 



in this work, the condition of ^> Vg,v r in Equation ()4]) is generally satisfied inside the 
disk. If the disk is perfectly in a stationary state, which means the centrifugal force (~ pvT) 
balance with the gravitational forces (e.g., Equation ([6])), no GWs can be emitted. As 
will be explained later, the disk attains mass continuously due to mass-accretion whose 
specific angular momentum increases outward (e.g., Figure 1). This is the primary reason of 
generating non-zero matter GWs from axisymmetrically but dynamically rotating collapsar 
disks. In the following computations, we assume that the observer is located in the equatorial 
plane (6 = tt/2 in Equation (P2])), and also that the distance to the GW source is comparable 
to nearly GRB-associated core-collapse supernovae (R = 100 Mpc) unless stated otherwise. 



2.2.2. GWs from anisotropic neutrino emission 



To compute the gravitational wavef orms from anis o tropic neutrino r adiati on, we follow 
the formalism pioneeringly proposed by lEpsteinl (119781 ); iMiiler fc Jankal fjl997l ). In the case 
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of our 2D axisymmetric case, the only non-vanishing component is the plus mode for the 
equatorial observer, 



AG 



dt 



d& $( 



y . dl v (9',t') 



c A R J " J ~ y ~ ' dVt' 
where $(#') depends on the angle measured from the symmetry axis (9 ) 

$(#') = 7rsinfl'(-i + 2| cos 6»'|). 



(9) 



(10) 



As given in Figure 1 of iKotake et al.l (120071 ). this function has positive values in the north 
polar cap for < 9' < 60° and in the south polar cap for 120° < 9' < 180°, but becomes 
negative values between 60° < 9' < 120°. To determine the anisotropy in neutrino emission 
(i.e., dl v /dVL in Equation Q), we perform a ray-tracing analys is, whic h was first proposed 



to be applicable in the Newton ian gravity (Kotake et al. 2009bl a. 2011 ) and later improved 



to be treatable in SR and GR ( Harikae et al.ll2010bl ; 



Harikae et al.ll2010al ). 



Applying the formalism in iLindquistl (119661 ). the Boltzmann equation for the neutrino 
occupation probability f u (e u ,Q) for a given neutrino energy (e u ) along a specified direction 
of Q can be expressed as, 



df u (e u ,Q) 
dX 



n[Q e (l - fu) - k/ u ] = n[Q e - K*f u ], 



'IV 



where n{x) is the proper number density of the external medium with which neutrinos 
interact and thus measured in its own local rest frame, e v is the neutrino energy measured 
in the local proper frame, A denotes an affine parameter along the geodesies, Q e and n 
represents neutrino emiss i vity an d absorptivity, and (1 — /) represents the Pauli blocking 



term (e.g., lHarikae et al.l (j2010al ) for more detailed information to derive the equation). 



Here we consider only v e and v e for simplicity, and the energy and momentum transfer via 
neutrino scattering are neglected, which is not only difficult to be treated by the ray-tracing 
technique but also a major undertaking in the radiative transport problem in general. In 
the final expression of Equation (fTTl) . n* is defined to represent the effective absorptivity. As 
for the opacity sources of neutrinos, electron capture on proton and nucle i, positron capture 



on neutron, neutrino scattering with nuc leon and nuclei, are included (IFuller et al.l 11985 
Takahashi et al.l 1 19781 ; iBruenn et al.ll2010l ). Here k* is estimated as k* = S[n targct • ^"(e^)] 



with n targGt , a{e u ) being the target number density of each reaction and the corresponding 
cross sections, respectively. 



According to IZinkl (120081 ). the formal solution of Equation (fTTT) can be given as 



/"Aout 

f v (e,Q)= / n(X")Q e (X",e u )exp 



n{X')K*{X')dX' 



A" 



dX" 



(12) 
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which is referred to as the rendering equation of the radiation transport problem. We perform 
a line integral along the geodesies from every point on the surface of neutrinospheres (A ) 
from which neutrinos can escape freely, up to the outer-most boundary of the computational 
domain (A out ). On the neutrinospheres, the neutrino distribution function is assumed to 
take a Fermi-Dirac type f FD [= l/(e u /k^,T + 1)] with a vanishing chemical potential For 
the neutrino energy bins (e„), we use 16 logarithmically spaced energy bins reaching from 3 
to 300 MeV. 

The path integration in Equation (I12p is done explicitly along the geodesies. In do- 
ing so, we determine each integration step by restricting the maximum change of neutrino 
opacity for all the neutrino energy-bins to be less than 10 %. By using an adaptive-mesh- 
refi nement approa c h also, our ray-trace code was proved to safely pass several test problems 



sec 



Harikae et al.l (l2010af ) for more details). For the ray-tracing calculation, 300(r) x 40(9) 
x 32(0) meshes points are cast. The radiation field typically changes more slowly than the 
hydrodynamics. We perform the ray-tracing calculation every 50 hydrodynamic step to re- 
duce the computational cost. Note that the obtained waveforms did not change significantly 
when the interval is varied as every 100 or 200 step. 

With f(e u ,Q) at the outer-most boundary, which is obtained by the above procedure, 
the emergent neutrino energy fluxes along the specified direction of Q can be estimated, 

dl v (n,e„) f e 2 v de v 

f(e u , Q) • (ce„) • (13) 



dfldS J J v ' v v) (2tt^c) 3 ' 

By summing up the energy fluxes with the weight of the area in the plane perpendicular to 
the rays (:dS), we can find dl u /dVL along the specified direction Q, 

tup) _ f duo, dS 



dtt J dfldS 

Repeating the above procedures, dl u (Q)/dQ can be estimated for all the directions, by which 
we can find the amplitudes of the GWs from neutrinos through Equation ()9]). 



3. Results 

First we pay attention to the properties of GWs as well as the hydrodynamic features 
in models without magnetic fields. 



5 &b is the Boltzmann constant 
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Fig. 2. — Gravitational waveforms from the sum of neutrinos and matter motions for all the 
non-magnetized models. Note that a collapsar is assumed to be located at the distance of 
100 Mpc. 



Table 1. Model Summary 



Model At(s) |/Wx|(lCr 22 ) |h„|(10- 22 ) E^(M @ c 2 ) 



J0.6 


3.41 


0.47 


0.14 


2.09 x 10 


J0.7 


5.22 


1.15 


0.37 


5.12 x 10 


J0.8 


11.0 


2.47 


1.80 


3.91 x 10 


J0.9 


11.0 


3.25 


1.22 


5.72 x 10 


J1.0 


11.0 


2.86 


1.07 


4.45 x 10 


Jl.l 


11.0 


2.54 


0.98 


3.56 x 10 


J1.2 


11.0 


2.23 


0.86 


3.26 x 10 


J0.8B10 


2.97 


0.18 


0.037 


1.73 x 10 


J0.8B11 


0.53 


0.028 


0.0065 


9.00 x 10 



Note. — At represents the simulation time. |/x ma x| represents the 
total GW amplitudes at maximum during the simulation time, while 
h v is the neutrino-originated GW at the end of the simulations. 
is the total radiated energy in the form of the GWs in unit of Mqc 2 . 
The distance to the source is assumed to be 100 Mpc. 
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Fig. 3. — Similar to Figure [2] but for the gravitational waveforms only from anisotropic 
neutrino emission (left panel, e.g., Equation flS])) and only from matter motions (right panel, 
e.g., Equation (jEJ)). Note that only the selected models are drawn in this figure not to make 
it messy. 



Figure [2] shows the total GW amplitudes for all the computed models without magnetic 
fields (see also table 1 for a model summary). As seen, the total GW amplitudes generally 
show a positively growing feature with time. In addition, a sudden disappearance in the 
signals can be seen for some slowly rotating models (models JO. 6 (red line) and JO. 7 (green 
line)). Figure |3] shows the gravitational waveforms contributed only from anisotropic neu- 
trino emission (left panel) and only from matter motions (right panel), respectively. The 
positively growing trend is shown to come from the neutrino contribution (left), which is 
much larger than the matter contribution (right). To understand these properties, we first 
briefly summarize a hydrodynamic evolution of our 2D non-magnetized models in section 
13.11 And then in section I3.2[ we move on to analyze the reason of the positive growth in 
the neutrino GWs by performing the ray-tracing analysis. After that, we analyze the matter 
GWs paying particular attention to the magnetic effects in section l3~3l We then discuss their 
detectability in section 13.41 



3.1. Hydrodynamic Features 

To capture hydrodynamic features in our models, Figure H] shows the evolution of neu- 
trino luminosities for some selected models. In t ~ 1.0 s after we start our simulations (t = 
s), all the models experience rapid infall and the subsequent shock formation in the center, 



-14- 




Fig. 4. — Time evolution of neutrino luminosities (the sum of all the neutrino species) for 
models JO. 6, JO. 75, JO. 8, JO. 85, J1.0, and, J1.2, respectively. Note that the time is measured 
from the epoch when the simulations are started (i.e., the onset of gravitational collapse). 

which leads to a drastic increase and the subsequent decrease in the neutrino luminosities^. 

Here we take model JO. 8 as a reference, because the precollapse angular momentum 
is adjusted to be closest to the original progenitor. The top left panel in Figure [5] show 
a snapshot at t — 0.42 s near the shock formation. Note that the accretion mass in this 
epoch is typically greater than 2~ 3M Q in the center. The maximum mass of the neutron 



star of the Shen equati o n of state is in th e same mass range (IShen et all (119981 ). see also 
O'Connor k Ottl (hoilh : iKiuchi & Kotake! J2008h . And recent full GR simulations show 
that the mass o f the central object just before collapsing to a BH is typically < 2.3M 



(e.g., lOtt et all (120111 )) with its typical radius of several km (as inferred from their Figure 
3). These evidences might support our very crude assumption of the prompt BH formation 
that is modeled by setting the BH initially in the center, although such assumption can be 
only tested by full GR simulations using the same progenitor model. Later on, the density 
configuration (compare the left-half in each panel) deforms to be more oblate with time due 
to accretion of material with higher angular momentum outside (e.g., Figured]). As will be 
explained in the later section, the luminous accretion disk is the primary source of the GWs 



9 Note that the shock formation is not because the central density exceeds the nuclear density, but because 
the matter pressure becomes so high due to compression that it pushes back the ram pressure of accreting 
material in the vicinity of the inner boundary. I n con t rast to the state-of-the-art simulations of BH-forming 
core-collapse supernovae (e.g.. ISumivoshi et al.l (|2006f ) : iFischer et al.l (|2009l ); lOtt et al.l (|201ll) ). the inability 
of capturing dynamics correctly especially before the BH formation, is one of major drawbacks in collapsar 
simulations in general. 
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Fig. 5. — Three snapshots characterising the hydrodynamic evolution for model JO. 8. The 
pair panels (top left) are for t = 0.42 s, and the top right, bottom panels are for t = 1.89 
and t = 9.00 s, respectively. In the pair panels, the left half in the top panel shows the 
logarithmic contour of density (in g cm" 3 ), while the right half is for entropy per nucleon in 
unit of the Boltzmann const ant (/cb)- The velocity fields are drawn by the white arrows, and 
the length is normalized by the scale shown in top right edge of the box (i.e. 10 10 cm s _1 ). 
The central black circle represents the inner boundary of our computations. 
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from anisotropic neutrino emission. 

The top right panel in Figure shows a snapshot at t — 1.89 s. The bluish regions near 
along the spin axis of the accretion disk (left-half, density) correspond to the so-called polar 
funnel regions. The entropy becomes highest near the surface of the accretion disk due to the 
shock heating when the accreting material hits the wall of the disk (right panel). Comparing 
Figure [3] (left panel) to Figure HJ the GWs from neutrinos deviate from zero typically later 
than t > 2 — 3 s, when the (total) neutrino luminosities become as high as ~ 10 52 erg/s. 
Later on, the neutrino luminosities show a gradual increase with time, reflecting the increase 
in the mass accretion to the newly formed accretion disk. 
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Fig. 6. — Same as Figure [5] but for model JO. 6. The top left, top right, and bottom panel is 
for t = 3.09, 3.42, 3.43 s, respectively. For this model, the accretion disk is finally absorbed 
into the central BH due to its small initial angular momentum (bottom panel), leading to a 
sudden decrease in the neutrino luminosity (e.g., Figure @J. 



In Figure HI the neutrino luminosity for models JO. 6 (red line) and JO. 7 (green line) is 
shown to steeply decrease at ~ 3.3 s and 7.5 s, respectively. This is because the accretion disk 
is swallowed to the center (bottom panel in Figure [6]) mainly because the pressure support 
in the accretion disk is reduced by the neutrino cooling. The disappearance of the accretion 
disks is also the reason of the sudden decrease in the GW signals observed both in the neutrino 
and matter sectors (e.g., Figures [2] and [3]) . Before the disappearance, the accretion disk is 
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Fig. 7. — Same as Figure |3] but for model Jl.l. The top left, top right, and bottom panel 
is for t = 0.5, 2.4, 9.0 s, respectively. Due to a more larger initial angular momentum than 
that for Figure El the accretion disk is deformed to be more oblate. 
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observed to show a rapid expansion and contraction (as indicated by the top panels in Figure 
[6]). At the same time, the mass flux on the Lagrange point around the disk changes v i olently , 



which may be r elated to the so-called runaway instability (e.g., lAbramowicz et al.l ( 119831 ); 



Font fc Dai end (120021 ) ) . Except for these slowing rotating models, the neutrino luminosities 
gradually settle to be nearly constant typically later than t ~ 6 s (FigureSJ). The saturation of 
the neutrino luminosity is because the d i sk is a lready deleptonized and the neutrino emission 
there is suppressed (e.g., lHarikae et al.l (120091 ) for more details). The neutrino luminosities 



at this epoch become as high a s 1Q 52 ~ 53 erg s 1 , whi c h tou c hes the level of ~ 10% of t he 



accretion luminosities (see also IChen fc Beloborodovl ( 120071 ); ISekiguchi fc Shibatal ( 1201 ll )). 
Apparently the accretion disk is neutrino-cooling dominated. 

From Figure HI the neutrino luminosity of model JO. 8 is shown to be highest (e.g., 
blue line). Models with higher initial angular momentum have more extended disks due to 
larger centrifugal forces (compare the bottom panels in Figure [5] and [7]) , leading to lower 
density and temperature in the disks. This is the reason that the neutrino luminosities for 
models J1.0 (light blue line), J1.2 (yellow line in Figure H]) are lower in this order. Reflect- 
ing this, the GWs from neutrinos become highest for model JO. 8 (left panel of Figure [31 
see also \h u \ in Table 1). Regarding the matter GWs, it should be noted that their max- 
imum absolute amplitudes are obtained not for the most rapidly rotating model (model 
J1.2), but for mode JO. 8 (right panel of Figure [3]). This situation may be akin t o the 



GW si gnal s emitt ed near core-bounce in core-collapse supernovae (see iKotake et al.l ( 12006 



20111 ); lOttl (120091 ) for recent reviews). Too much initial angular momentum works to sup- 
press the matter compression leading to the smaller mass-quadrupole moment. As a re- 
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3.2. GWs from Anisotropic Neutrino Emission 

In this section, we move on to look more into detail the properties of gravitational 
waveforms mentioned in the previous section. By taking model JO. 8 as a reference, we first 
focus on the neutrino GWs. 

Figure [8] shows the local neutrino energy fluxes (:dl u / (dQdS), Equation ( 113]) ) at t = 9.0 
s (e.g., bottom panel of Figure [5]) seen from polar (left) or equatorial direction (right), 
respectively. Note that the polar direction is taken to be parallel to the spin axis of the 
accretion disk. The top panel is for the Minkowski geometry, which corresponds to a purely 
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Newtonian case. Reflecting the shape of the accretion disk in axisymmetry, the contours of 
the neutrino flux are deformed to be oblate when seen from the equator (top right), while the 
neutrino flux seen from the polar direction (top left) looks like a superposition of circles in 
a concentric fashion, the center of which corresponds to the spin axis. Note that the central 
small cavity in the top left panel corresponds to the inner boundary of our SR simulation. 

The middle pa nels are calcula t ed for the Minkowski geometry including special relativis- 



tic corrections (see lHarikae et al.l (j2010al ) for more detail). Comparing the middle right to 
the top right panel, a clear difference is the left-right asymmetry in the middle right panel. 
This is due to the special relativistic beaming effects. Since the rotational velocity of the 
accretion disk is perpendicular to the polar direction, the special relativistic beaming effect 
suppresses the neutrino emission toward the polar region. As a result, the neutrino fluxes are 
made dark (right-hand side in the middle right panel) because the direction of the rotating 
material is opposite to the direction to the observer. To study the GR effects on the neutrino 
luminosity, we place 4M© in the central region, which mimics the event horizon of the BH. 
The effects of the BH spin are examined by setting the Kerr parameter by hand as a = and 
a = 0.999 for the Schwarzschild and the extreme Kerr geometry, respectively. The bottom 
panels of Figure M are for the extreme Kerr geometry, which however looks very similar to 
the middle panels. 

For a more detailed comparison, the left panel of Figure [9] shows the neutrino luminosity 
per solid angle dL v /dVL for the Minkowski (indicated by "SR"), Schwarzschild ("GR(a=0))", 
and extreme Kerr geometry ("GR(a=0.99))" , respectively. The most important message in 
this panel is that in every case, the neutrino luminosity seen from the direction near parallel 
to the spin axis (8 = 0) becomes higher than the one seen from the equatorial direction 
{9 = it/2). This is because the cross section of the pan-cake like accretion disk seen from the 
spin axis becomes larger compared to the one seen from the horizontal direction (Figure [S]). 
Remembering again that $(#') in Equation (1101) is positive near the north and south polar 
caps, the dominance of the polar neutrino luminosities make the neutrino GWs positive in 
the polar cap regions (right panel in Figure [9]). This is the reason that the neutrino GWs 
increase monotonically with time (e.g., left panel of Figure [3]). Here it is worth mentioning 
that from neutrino luminosity (L u in Figure H]), their typical duration (At u ), anisotropy in 
the neutrino radiation (a^, read from the left panel of Figure the GW amplitudes from 
neutrinos can be estimated from Equation as 

K ~ 10 " 22 (oi) (l053 erg s-0 (loi) (lOOMpc) ' (15) 
which is in good agreement with the numerical results (Figure [2]). 

Comparing to the SR with GR case in the left panel of Figure HU the neutrino luminosity 
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Fig. 8. — The neutrino energy fluxes of dl u /(dQdS) (Equation flTSJ)) seen from polar (left) 
or equatorial direction (right panels), respectively for model JO. 8 at t — 9.1 s. The top and 
middle panels are calculated for the Minkowski geometry without or with special relativistic 
corrections, and the bottom panels are for the extreme Kerr geometry (a = 0.999, see text 
for more details). Note that the Z axis in the right panels conincides with the spin direction 
of the accretion disk. 




Fig. 9. — The neutrino luminosity per angle (left panel, dl v /dVL (Equation (15)) and neutrino 
GWs per solid angle (right panel) for model JO. 8 at t =9.0 s in the case of the Minkowski (in- 
dicated by "SR"), Schwarzschild ("GR(a=0))", and extreme Kerr geometry ("GR(a=0.99))'\ 
respectively. 

becomes smaller by ~ 20 — 30% in the GR case due to the GR redshift and also due to the 
bending effects. Comparing the green line with the blue line, the frame- dragging effect due 
to the BH spin barely affects the emergent neutrino luminosity. This is probably because 
the GR effect on the neutrino lu minosity becom e s impo rtant only in the central regions very 



close to the BH (e.g., Figure 15 lHarikae et al.l (]2010al )). Our results indicate that a ray- 



tracing calculation in the Schwarzschild geometry is at least needed to accurately estimate 
the neutrino GWs in the collapsar's environment. 



3.3. GWs from Matter Motions and Magnetic Fields 

Now we shortly analyze the matter GWs in model JO. 8. As summarized in section 
12.21 the matter GWs are estimated by the sum of the hydrodynamic and gravity parts in 
Equations (jU [6]) for non- magnetized models. Among the kinetic energy terms (oc pViVj) in 
Equation (TJJ, the largest contribution comes from the rotational energy (i.e., —p if W 2 v^), 
which is shown in Figure [10] (green line). In contrast to the negative contribution by this 
term, the gravity part (blue line labeled by " Gravity" in Figure fTUl) is shown to make a 
positive contribution. If a star rotates perfectly stationary, the centrifugal forces balance 
with the gravitational force, leading to no GWs. In the collapsar disk studied in this work, 
the disk attains mass continuously due to mass-accretion with increasing its specific angular 
momentum outward. This is the reason why the disk is not perfectly stationary, leading to 
a non-zero GW emission from matter motions (see pink line labeled by "Matter" in Figure 
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Fig. 10. — The gravitational waveform for model JO. 8. "Total" denotes the total amplitudes 
(the sum of /i( ma tter) in Equation flH]) and h v in Equation OH])), while "Matter", "Rotational 
Energy", "Gravity", and "Neutrino" represents the contribution from matter (/t( ma tter) in 
Equation flH])), rotational energy (see text for the definition), gravity (fygrav) in Equation 
OH])), and neutrino emission (h u in Equation ©), respectively. 
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ITUl) . However, their GW amplitudes are much smaller compared to those from anisotropic 
neutrino emission (compare the pink line with the light blue line (neutrino GWs) in Figure 

HDD- 

The left panel of Figure [11] shows a normalized contribution of each term in Afo at 
i = 9 s for model JO. 8, which is estimated by the volume integral of Afo within a given 
sphere enclosed by certain radius. The region in a radius between 50 km to 160 km is shown 
to contribute to produce GWs, which corresponds to a high density region in the accretion 
disk (see the bottom panel in Figure E]). It can be also shown that the radial gradient of 
the gravitational potential (green line indicated by "Gravity (r)" in the plot) closely cancels 
with the contribution from the centrifugal force (red line) in the disk. As shown in the right 
panel of Figure HU the enclosed mass there becomes as big as ~ 4M as a result of the 
hyperaccreting activity lasting ~ 9 s till then, and the typical rotational period there is the 
order of milliseconds (red line). 



Putting these numbers to the standard GW stress formula (e.g., IShapiro &: Teukolsky 



( 119831 )). an upper bound of the matter GW amplitudes may be estimated as 
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where D is the distance to the source, Ijj is the second time derivative of the quadrupole 
moment of .7^, M, R, and T represents the typical mass and radius of the accretion disk and 
the timescale which we may take as the rotational period, respectively, and e is the degree 
of the nonsphericity, which we take optimistically as 10 % in the case of the accretion disk. 
This estimate, roughly consistent with the numerical results (e.g., right panel of Figure 3), 
also shows that anisotropic neutrino emission, producing the GW amplitude on the order 
of w 10 -22 for a source of 100 Mpc, is the primary source in the long term evolution of 
collapsars (~ 10 s) (see Equation ffT5l) ). 

Now we discuss the effects of magnetic fields on the waveforms by taking model J0.8B11 
that has a strongest initial magnetic field in our models. In Figure [T2| the pink line represents 
the GWs from magnetic fields (Equation Q). But, first of all, let us shortly comment on 
the burst-like feature of the green line (contributed from matter motions) at t ~ 0.24 s 
in the figure. This may look similar to the burst of GWs emitted at the moment of t he 
BH formation possibly followed by a damped sinusoid al oscillation (ISeidel dl990l 119911)) 



but this can be only captured in full GR simu lations f lBaiotti et al.l 12007c lOtt et al 



2011 



Sekiguchi fc Shibatal 12011c iKuroda et al.ll2012l ). As already seen in the luminosity plot of 
Figure HI the above burst simply corresponds to the shock formation at the center in our SR 
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Fig. 11. — Left panel shows a normalized contribution of each term in A^q (e.g., Equation 
([3]) as similar to Figure [TO]) as a function of radius for model JO. 8 at t — 9 s. "Gravity 
{r/9)" indicates contributions from radial and lateral derivative of the gravitational potential 
(corresponding to the last two terms in Equation (6)). Right panel shows the radial profiles 
of the rotational period (red line) and the enclosed mass normalized by M Q (green line), 
respectively. 



models. Back to the main point, we analyze the increasing trend of the GWs from magnetic 
fields in the following. 

The left panel of Figure [TBI depicts a snapshot of density (left-half), entropy (right-top) 
and plasma j3 (right-bottom, the ratio of magnetic to the matter pressure) at the final simu- 
lation time (t = 541 ms) for model J0.8B11. The high entropy regions in the slightly off-axis 
region (seen as reddish in the top right panel) correspond to the magnetohydrodynamically- 
driven outflows that are pushed outwards by the twisted toroidal magnetic fields (also seen 
as reddish in the bottom right panel, indicated by "Toroidal"). 

The right panel of Figure [13] shows contributions to the total GW amplitudes (Equation 
(9)), in which the left-hand-side panels are for the sum of the hydrodynamic and gravitational 
part (indicated by "Matter"), namely log ^± ^4f 2 ( hyd ) + ^lo(grav) J (^ top(+)/bottom(— ) 
(Equations (]4l6]) ). and the right-hand-side panels are for the magnetic part, namely log 
[^±Afo (mag) J] (right top(+)/bottom(— )) (e.g., Equation (J7D). By comparing the top two 
panels, it can be seen that the positive contribution is overlapped with the regions where 
the MHD outflows exist. The major positive contribution is from the kinetic term of the 
MHD outflows with large radial velocities (e.g., +p*W 2 v r 2 in Equation (j3|). The magnetic 
part also contributes to the positive trend (see top right-half in the right panel (labeled by 
mag(+))). This comes from the toroidal magnetic fields (e.g., +b < f > 2 in Equation (|7j)), which 
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dominantly contribute to drive MHD explosions. 



Unfortunately, our MHD code becomes numerically unstable when the strong MHD 
jets propagate to a stellar mantle with decreasing density, which prevents us from studying 
the resulting GWs in more long run. The neutrino GWs are much smaller than the other 
GW sources (e.g., Figure [12]) ) simply due to the shorter simulation time. We expect that 
the increa sing trends by the magneti c fields maintain as the MHD shocks propagate further 



out (e.g., iTakiwaki fc Kotakd ( 120111 )). To confirm it, we need to implement a numerical 



technique special ly developed to solve the force-free fields, which is major undertaking (e.g., 
McKinnevI fco06h \ 




Fig. 12. — Similar to Figure |2] but for model J0.8B11. The pink line represents the GWs 
from magnetic fields (/i(ma g ) in Equation flB])). 



3.4. Detectability 

Finally, Figure [Ml depicts the GW spectra for models JO. 8 (left) and J1.0 (right). The 
GW spectra in the frequency domain between 1 to 100 Hz becomes slightly larger for model 
JO. 8 than for model J1.0. This reflects a more efficient release of the gravitational bind- 
ing energy for the moderately rotating case (model JO. 8) as already mentioned in sections 
3.2 and 3.4. For the cosmological distance scale of GRBs (~ 100 Mpc), these low fre- 
quency GW signals are unfortunately very hard to detect even by the advanced detectors 
(like the advanced LIGO or KAGRA/LCGT) whose sensitivity is severely limited by the 
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Fig. 13. — Hydrodynamic configuration for model J0.8B11 at t =541 ms. In the left panel, 
the logarithmic density (in g cm -3 , left-half), entropy^s/baryon, right top), and the ratio 
of magnetic to the matter pressure (right bottom) are shown. Right panel shows the sum of 
the hydrodynamic and gravitational parts (indicated by "Matter" in the left-hand side) and 
the magnetic part (indicated by "Mag" in the right-hand side), respectively. The top and 
bottom panels represent the positive and negative contribution (indicated by (+) or (-)) to 
A|q, respectively (see text for more detail) 
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Fig. 14. — Spectral distributions of GWs from matter motions for models J0.8 ( left) and J1.0 



(right) with the expected detection limits of advanc ed LIGO flWeinsteinl 120021), Large-scal e 



Cryogenic Gravitational wave T elescope (LCGT) (IKuroda fc 



CGT Collaboration! |2010| ) 



and Fabry- Perot type DECIGO (lKawamurall2006l ; iKudoh et al.l 120061 1 . The distance to the 
so urce is assumed to be lOOMpc. /z c h ar is the characteristic gravitational wave strain defined 



m 



Flanagan fc Hughes! (119981 ). Note that the frequency domain larger than 0.1 Hz is only 



plotted, since the typical timescale that our simulations covered is 10 s. 
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seismic noises (lAndo fc the TAMA Collaboration! bOQa lAbbott et al.ll2005l : IWeinsteinl 12002 



Kuroda fc LCGT Collaborationll2010l ). A good news is that these signals could be detectable 



by a recently prop osed future space-based interferometer like Fabry-Perot type DECIGO 
( jKawamural ( 120061 ). black line in Figure lT4"j) . Two low- luminosi ty LGRBs were already ob- 
served at distances of ~ 40 Mpc (98042f UGal£una et all (jl998h ) and ~ 130 Mpc (060218, 



Ferrero et al.l (120061 ) ). And their local rate, being much higher than that of normal bursts, 



is expected to be as large ~ -l.Dmn yr" 1 with ffinn re presenting the distance normalized 
by 100 Mpc (see discussions in lCorsi &: Meszaroj (120091 )). Our results suggest that the GW 
astronomy of collapsars could become reality by the DECIGO-class GW detectors, hopefully 
near in the future. 
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4. Summary and Discussion 



By performing axisymmetric SRMHD simulations, we investigated possible signatures 
of GWs in the context of the collapsar model of LGRBs. By cutting out the central BH, 
we focused on the GWs generated by asphericities in neutrino emission and matter motions 
in the vicinity of the hyperaccreting disks. Nine models were computed by changing initial 
angular momenta and magnetic fields parametrically in the precollapse core of a 35M Q 
progenitor star. As for the microphysics, a realistic equation of state was employed and 
the neutrino cooling was taken into account via a multiflavor neutrino leakage scheme. To 
accurately estimate GWs from neutrinos, we performed a ray-tracing analysis in GR by a 
post-processing procedure. We studied also the effects of magnetic fields on the gravitational 
waveforms by employing a stress formula that includes contributions both from magnetic 
fields and special relativistic corrections. We found that the GW amplitudes from anisotropic 
neutrino emission shows a monotonic increase with time, whose amplitudes are much larger 
than those from matter motions of the accreting material. We showed that the increasing 
trend of the neutrino GWs stems from the excess of neutrino emission in the direction near 
parallel to the spin axis illuminated from the hyperaccreting disks. We pointed out that a 
recently proposed future space interferometer like Fabry-Perot type DECIGO would permit 
the detection of these signals within rs 100 Mpc. 

Here it should be noted that our 2D simulations cannot capture any non-axisymmetric 
instabilities so far proposed to provide a strong GW emis si on in the semi-analytic a l (e.g ., 
van Puttenl f200lh: E)avies et all (12002 ) ; iFrver et all J2002h ; iKobavashi fc Meszarosl (120031): 
Piro fc PfahlJ2007h:ICorsi & Meszarosl J2009h ) and in full GR simulations (e.g Jshibata fc Sekiguchi 



( 20051 ); Manca et al. ( 2007 )). Therefore the present results might be regarded to give a lower 
limit for the possible GW emission in collapsars. To go up the ladders beyond the 2D sim- 
ulations is very numerical challenging, however, we need to handle it not only to test the 
outcomes of the proposed ideas about the non-axisymmetric instabilities but also to obtain 
more accurate waveforms from collapsars. 

When studying the formation of BHs and the associated GW signals, the use of a 
pseudo-newtonian potential can lead to significant errors. This is the reason why we had to 
limit our discussion only to the asphericities and the resulting GW emission in the vicinity of 
the accretion disk which is far away from the central object. In studying the dynamics from 
core-collapse to neutron star, a conformally -flat- condition (CFC) approximation has been 
often employed to solve the GR equations (IDimmelmeier et al.ll2002l . 120071 120081 ) . And it 
has b een tested that such a treat ment is very good in capturing the results of full-GR results 
(e.g., IShibata fc Sekiguchi! ( 120031 1) for a wide variety of supernova progenitors with neutron 
star formations. But it may need a further investigation that the approximation is still valid 
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in applying to collapsars which is a highly aspherical disk-BH system with very dilu t e pola r 



funnel reg i ons a l ong the rotational ax i s. Fu ll GR simulations (e. g., iBaiotti et al.l ( 120071 ); 



Ott et all (1201 If ): ISekiguchi fe Shibatal (1201 If ): Kuroda et al.l ( 120121 ) and references therein) 
are indeed one of the most important topics pointing to the final frontiers of stellar core- 
collapse simulations, however, it is generally computationally too expensive at present to fol- 
low the late-time collapsar evolution up to ~10 s, which is a typical duration of long-duration 
GRBs. Moreover the inclusion of microphysics such as neutrino heating is a major undertak- 
ing for the full GR simulat i ons. T h e Cowling approximation ( or the fixed metric approach 



McKinney fc Narayanl (]2007bf ); iKomissarov fc McKinneyl (120071 ) ; lBarkov fc Komissarov 



( 120081 )) has been often used in collapsar simulations so far. But it still remains a non-trivial 
issue how to treat the self-gravity of the accretion disk. Finally, the pseudo-newtonian ap- 
proach we take in this work cannot unambiguously capture accurate properties of the flows 
in the vicinity of the BH as well as the associated GW signals. Sacrificing the central regions, 
such simplified method would be currently only a possible way to follow a long-term evolution 
(especially for the disk evolution) with including an appropriate treatment of microphysics. 
Needless to say, these four approaches (full GR, CFC, Cowling, and post-newtonian) may be 
regarded to be complimentary useful to study the different epochs in the collapsar evolution. 
For bridging the gaps between them, i t may be a good idea to employ the end results of fully 
GR simulations (e.g., lOtt et al.l (120111 )) in our simulation, which we consider to be the most 



urgent task to investigate as a sequel of this study. 

In a data analysis to extract the true GW signals from the confusing detector noises, it 
should be of primary importance to take a coincident analysis with the conventional electro- 
magnetic observation as well as neutrinos. What could be the photon and neutrino signatures 
in our collapsar models ? To answer this question, we need to perform a long-term simulation 
that bridges continuously the phase 2 and 3 as mentioned in the introduction. For generating 
neutrino-driven or MHD-driven outflows in numerical simulations of collapsars, the GRMHD 
simulations including the e ffects of the neutrino h eating are needed, t h e formulation of whi ch 
is in a steady progress (e.g.. lShibata et al.l ( 1201 ll ) ; iMiiller et al.l ( l2010l ); lKuroda et al.l ( 120121 )). 
Assisted by a growing computational power, these advanced simulations would be hopefully 
practicable by utilizing the next-generation supercomputers. These updates should bring 
forward not only our understanding the dynamics of the collapsar en g ines, but also the 
theoretical predictions of observable multi- messengers (e . g.. lAndo et al.l (120121 ) for a recent 



revie w 10 !), including neutrin o emission (e.g.. lAbbasi et al.l ( 1201 ll )) and nucleosynthetic yields 



(e.g., iFujimoto et al.l (120071 ) and references therein). The neutrino signals emitted from our 
collapsar models (e.g., from Figure H] with its typical emergent neutrino energy of ~ 15 — 20 



see, iKotake et al.l (|20 12T ) for a review about multi-messenger perspectives on core-collapse supernovae. 
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MeV) are visible up to the Local Group (~ 1 Mpc, Kawagoe et al. in preparation) by fu- 
ture megaton-class detect ors (e.g., Hyper-K amiokande, Memphys, and LBNE), large-scale 



scintillators (e.g ., HALO (lEngel et al.ll2003l ). and by liquid- Argon detectors (GLACIER, see 
Scholbergl ( 120101 ) for collective references therein). While this work might raise many more 
questions than it can answer, it definitely makes clear that our understanding of the GWs 
in collapsars is still in its infancy and that collapsars and the BH-forming supernovae are 
"gold mine" in which a number of unsettled and fascinating research themes are hidden. We 
hope that our exploratory results, at least, give momentum to theorists to make the GW 
prediction based on a more sophisticated numerical modeling of collapsars. 
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